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We model charged colloidal suspensions using an analogy with foams. We study the solid-solid 
phase transitions of these systems as a function of particle volume-fraction and ionic strength. The 
screened Coulomb interaction is replaced by an interaction between walls of the Voronoi cells around 
each particle. We fit the surface charge to reproduce the phase diagram for the charged suspension 
studied by Sirota et al. [Phys. Rev. Lett. 62, 1524 (1989)]. With this fit parameter we are able to 
calculate the elastic moduli of the system and find good agreement with the available data. 
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Charged colloids have been the subject of intense 
study both experimentally and theoretically. In the 
laboratory, optical techniques can readily probe these 
systems and they are easily manipulated via chemical 
means ^, ^, Q]. Their rich chemistry leads to many 
industrial applications ranging from uses in coating ma- 
terials, in ceramic precursors, and in designing and man- 
ufacturing biological macromolecules. These systems in- 
teract through a short-range screened Coulomb poten- 
tial and thus epitomize one of the few classes of interac- 
tions in soft systems. In this Letter, we replace the exact 
screened-Coulomb interaction between colloidal spheres 
with an effective interaction proportional to the surface 
area of the Wigner-Seitz or Voronoi cells which contain 
each sphere. By balancing the entropic interaction which 
favors the close-packed face-centered cubic (FCC) lattice 
with the surface interaction which, as we will argue, fa- 
vors the body-centered cubic (BCC) lattice, we will find 
the phase coexistence line as a function of volume fraction 
and salt concentration. We compare our results to exper- 
iment to find the surface charge on the spheres. Using 
this fitted value, we calculated the shear moduli within 
our framework and found close agreement with experi- 
ment. Our model leads to quantitatively different results 
than those found in molecular dynamics (MD) simula- 
tions 1^ of particles interacting via Yukawa interactions. 

Several experimental studies ||l|, ||, || provide a wealth 
of data on the stability of colloidal phases, in particu- 
lar, the fluid phase, the glass phase and two solid phases 
with either FCC or BCC lattices. These systems are 
aqueous suspensions of uniformly charged polystyrene 
spheres with a variable salt concentration which is a con- 
trol parameter for the degree of screening of the underly- 
ing Coulomb interaction. Because the lattice spacing is 
large in the ordered phases, they can be probed by scat- 
tering of visible light Q . Through scattering it is possible 
to determine those volume fractions at which the order- 
disorder transition occurs at low salt concentrations and 
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at which the FCC-BCC transition occurs. Though the 
MD simulations qualitatively corroborate these exper- 
imental findings, one is hard-pressed for similar quantita- 
tive agreement between existing experimental data and 
simulation results. The lack of exact analytic methods 
with which to compare experiment or simulation further 
complicates this problem. Thus theoretical models which 
distill the crucial aspects of the interactions and lead to 
testable predictions are necessary. 

Recently, a foam analogy was used to account for the 
many non-closed packed structures observed in colloidal 
and lyotropic systems 0| . In this analogy, the surface 
interaction replaced the brush-like interaction between 
the spheres and accounted for the entropic degrees of 
freedom associated with configurations of the soft coro- 
nas of alkyl chains attached to hard cores. Since the 
coronas prefer a packing arrangement that disfavors in- 
tcrdigitation and thus maximize the average separation 
of the cores, lattices with small interfacial area are fa- 
vored in such systems. For example, the BCC lattice 
of orthic tetrakaidecahedra has a smaller area than the 
FCC lattice of rhombic dodecahedra, and so the surface 
interaction favors the BCC lattice. 

Here we apply the same framework to understanding 
the experimental phase diagram |^ and MD simulation 
results for charged colloidal systems. In particular, 
our model couples the bulk free energy of the charged 
hard cores to the free energy of their screened Coulomb 
interaction via the effective surface area. The model con- 
tains a single free parameter: the surface potential asso- 
ciated with each colloidal particle. Using experiment as 
our guide, we fit the surface potential to reproduce the 
phase diagram. In order to test our theory, we use the 
same surface potential to calculate the shear and bulk 
moduli and compare these with experiment and simula- 
tion. We have found good agreement with the available 
data and, in doing so, demonstrate the plausibility of un- 
derstanding the various experimental results of charged 
systems in terms of a geometrical principle already es- 
tablished for another class of colloids. 

An apparent dramatic feature of our results is the neg- 
ative values of the bulk modulus for small A — na, where 
K is the inverse Debye screening length and a is the av- 
erage interparticle spacing. This is an artifact of our 
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model which can be simply understood: for A = 1 the 
Debye screening length is exactly equal to the interparti- 
cle spacing. As a result, we would expect next-to-nearest- 
neighbor interactions to play an important role in the sta- 
bility and energetics of the charged colloidal crystal and 
therefore would not trust our model at small A. However, 
for larger values of A the bulk modulus become positive 
- for the data we compare with, A ~ 4.64, well within 
the regime where next-to-nearest-neighbor interactions 
are small. Our model also yields qualitatively different 
behavior for the relative magnitudes of the elastic con- 
stants in comparison to MD simulations ||] . We attribute 
this difference to the finite size of the colloids - in the MD 
simulation they were treated as point particles and thus 
did not suffer entropic losses from excluded volume. We 
also note that the shear modulus becomes negative at 
large values of A. This signals the instability of the BCC 
lattice and its transition to FCC: for large A the screened- 
Coulomb potential becomes unimportant and the system 
should behave as if the interactions were purely hard- 
core. This is in line with the recent experimental de- 
termination 1^ that the FCC lattice is the equilibrium 
phase for samples with A > 10. Moreover, in that study 
it was argued that since Yukawa-like pair interactions 
could not account for both the elasticity of the crystal 
and the colloidal dynamics, three-body and higher inter- 
actions must play an important role. Since our approach 
is intrinsically many-body, our results for the elastic con- 
stants as a function of A are consistent with that conclu- 
sion - we find behavior that cannot be explained by MD 
simulation of particles interacting via Yukawa potentials. 
There have been alternative explanations of this discrep- 
ancy based on hydrodynamic screening jl^. We hope 
that our framework will allow for improved understand- 
ing of these issues. 

In a colloidal suspension at fixed density, the volume of 
the container is the sum of the volume of the hard spheres 
and the excess volume of the salt solution. Since this 
latter volume can be viewed as enveloping the individual 
spheres, we can imagine breaking the volume up into a 
lattice of Voronoi cells, each of which contains a colloidal 
particle. The excess volume can then be written as the 
product of the area of these dividing surfaces and their 
average thickness so that 

Ad — constant, (1) 

where A is total area of these bilayers and d is their av- 
erage thickness. Maximizing d to reduce the repulsive 
interaction amounts to minimizing A whence comes the 
interaction that favors BCC lattices (and the more exotic 
A15 lattices |,|lll). Since the crystal lattice is composed 
of identical colloidal particles in cases of interest, the min- 
imum area problem is equivalent to Kelvin's problem fl^ 
of partitioning a given space into cells of equal volumes 
having the smallest interfacial area. 

To compute the bulk free energy of the system, we fol- 
low the method developed in a previous study ^ by 
adopting the cellular free-volume theory. In this approx- 



imation, each particle is confined to a cage formed by its 
neighbors. The free volume available to each particle's 
center of mass is the volume of the Wigner-Seitz cell af- 
ter a layer of thickness (t/2 (where cr is the hard-core 
diameter of particles) is peeled off of its faces. Despite it 
mean-field nature, the free volume theory yields excellent 
quantitative agreement with available numerical simula- 
tions in the high-density limit and is roughly 1% off above 
the hard-sphere fluid-solid transition. The use of free vol- 
ume theory at the lower densities probed in the charged 
colloidal systems could be suspect. However as long as 
the shear elastic constants are non-zero, we expect that 
an "Einstein-crystal" description of the phonon modes 
should be adequate. Indeed, we will find that appropri- 
ate moduli are non-zero in the density and salt concen- 
tration regime we are studying [ p^ and therefore there 
are no "soft modes" which might strongly contribute to 
collective effects. The resulting bulk free energy for the 
FCC or BCC lattice is: 

Fb-ik = -ksTln(a- (rn^Va _ (2) 

where n — pa^ is the reduced density and a is the hard- 
core diameter of each colloidal particle. The coefficients 
Q,Fcc _ 25/2 g^j^jj Qf^'^'^ = 6.716 depend on the shape of 
the ceUs, whereas P"'"' = 2^1^ and = 2-2/331/2 are 

determined by their size. 

To obtain the electrostatic energy between the charged 
spheres, we model their interaction as that between two 
charged flat plates using the Debye-Hiickel approxima- 
tion |l4| for a salt solution with counterions of charge Ze. 
We replace each sphere with its Voronoi cell and then, in 
a Derjaguin-like approximation, replace the interaction 
with an interaction across parallel plates. The energy be- 
tween two plates can be found via an electromagnetic po- 
tential satisfying the linearized Poisson-Boltzmann equa- 
tion. For the experimental system in question we will 
show that this approximation is valid. Using linear su- 
perposition and the solution for a single charged plate, 
we have ||]: 

Fc = 64AMA:BTnbK"Hanh^ (i*'') "^^^ ' 

Taken together, Eqs. (||) and (||) encode the complete 
picture of the structural FCC-BCC transition in terms of 
the bulk countcrion number density nb, the dimcnsionless 
surface potential of colloids ^'s, the co lloid density n , 
and the Debye screening length kT^ — ekgT /2e'^Z'^n}j, 
itself a function of these control variables (where e is the 
dielectric constant). In order to relate the interaction 
energy to the reduced density n, we relate the spacing 
of the charged plates to the density via the constraint 
Eq. (^. The interparticle separation d depends on the 
particle density: 
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where = 7^cr^n~^/^, and 7^<^<^ = 5.345 and 

ryBCC _ 5.308 are dimensionless quantities characteriz- 
ine; the magnitude of area. Substitution of Eq. into 
Eq. (H) completes the expression for the interfacial free 
energy. Since the density and salt concentration are mea- 
surable parameters, the only free parameter is the surface 
potential ^s- It is therefore customary to treat as a 
parameter chosen for the best fit to the data |0. From 
the surface potential we will be able to determine the net 
surface charge q on the colloidal spheres. 

Focusing on the structural phase transition from BCC 
to FCC, we rely on the X-ray scattering study of charged 
polystyrene spheres by Sirota et al. That experi- 
ment studied the phase diagram of polystyrene spheres 
in a 0.9-methanol-O.l-water suspension. The thermody- 
namic behavior of this system was recorded as a function 
of volume fraction and salt concentration chci • In order 
to determine the surface potential, we choose an FCC- 
BCC coexistence point in the experimental data, and 
then equate the corresponding FCC and BCC free en- 
ergies to solve for Vfs- Through Gauss's law, we calculate 
the charge per colloid q — 2AM\/2efc^TVib sinh(ivl/g). 
Finally, using this value of the charge, we find the den- 
sity n at which F^^^ /kgT and F^^^ /ksT are equal for 
each salt concentration nb. Due to the limited data, there 
is only one experimental point at non-zero salt concen- 
tration that is in the coexistence region. We thus obtain 
the coexistence curve in ((/<-nb) phase diagram (Fig. |^). 
The relevant experimental data from Ref. [B| are included 
for comparison. 
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FIG. 1: Theoretical FCC-BCC coexistence curve as a func- 
tion of volume fraction <j) and electrolyte concentration (HCl) . 
The diamonds are the theoretical predictions, while the other 
points come from the data of Ref. Q . Solid squares are coex- 
istence points, open triangles are FCC, and open circles are 
BCC. 

Within our framework we find a surface potential of 
^'s = 0.2, and so we were justified in linearizing the 
Poisson-Boltzmann equation. However, this potential 
corresponds to a total charge of about 48e per colloidal 
particle. The value of the charge is smaller than the 



quoted experimental value of about 135e per sphere [H, 
though it is comparable in magnitude. The experimen- 
tal value was actually determined indirectly by measur- 
ing the shear modulus for the crystalline sample at 12% 
volume fraction | [l5[ |. Because of the uncertainty in the 
charge, we will make a direct comparison with experi- 
ment by calculating the shear modulus of the BCC lattice 
using our cellular framework. 

In general, the elastic energy of a cubic crystal has 
three elastic constants, Kii,Ki2, and ^^44 [16|: 



F, 



cubi( 



(5) 



For a polycrystalline sample, the bulk modulus K = 
(K\\ ^ 1K\2) /3 is an isotropic quantity. On the other 
hand, the shear modulus depends on the direction of 
the applied shear and ranges from the elongational shear 
/i = K\\ — K12/2 to the simple shear Ka^a- We determine 
the three elastic constants by calculating the free volume 
and the change in surface area for a deformed BCC lat- 
tice. The free volume was found via a simple numerical 
scheme while the change in surface area was found with 
the aid of Surface Evolver ||l7|. Equating the contin- 
uum energy change from Eq. ^) with the energy change 
within our framework, we find that Kn « 17.1 N/m^, 
K12 w -1.2 N/m3, and K^i w 10.3 N/m^ for a 12% 
volume-fraction (n = 0.23) BCC sample. Thus the shear 
modulus ranges from 10.3 to 17.7 N/m^. Again this is 
the same order of magnitude as the measured, isotropic 
value of 17 N/m^ which is a complicated combina- 
tion of the two shear moduli. We have calculated these 
moduli using an HCl concentration of 50 /imol, as this is 
in the middle of the reported BCC regime 

In order to compare with MD simulations we tab- 
ulate (Table ||) values of the bulk modulus and the three 
elastic constants as a function of A = Ka (in the exper- 
imental study, A w 4.64). However, unlike point parti- 
cles with Yukawa interactions, our model does not scale 
simply with this parameter. To vary A, we keep the col- 
loidal density fixed as well as the total charge on each 
colloid. Varying the salt concentration changes both k 
and the surface potential ^^s. With these new values we 
can recalculate the elastic constants. For completeness, 
we tabulate the values of the bulk modulus K and elastic 
constants of the BCC structure a,t (p — 0.12 for different 
values of A = Kfl (where a is the mean separation of 
spheres; at n = 0.23, a « 140 nm): 

While our framework assigns a surface charge that is 
smaller than the value suggested in the original analy- 
sis , it is in fact quite common for the effective charge 
to be reduced when compared with the net or titratable 
charge. To calculate the renormalized charge Z* and 
acid concentration n^Q^ , Alexander et al. employed 
a Wigner-Seitz approximation similar in principle to our 
way of determining the bulk free energies of the various 
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A 



K 



Kl2 



1 -26.1 / -7.59 17.5 / 10.7 -47.9 / -16.7 54.4 / 28.4 

2 -7.54 / 13.8 19.0 / 18.4 -21.1 / 11.6 29.1 / 15.6 

3 0.12 / 22.3 19.7 / 21.4 -9.67 / 22.7 18.3 / 10.6 

4 3.66 / 26.2 18.4 / 22.4 -3.71 / 28.0 12.7 / 7.95 

5 5.37 / 28.0 16.6 / 22.5 -0.22 / 30.7 9.40 / 6.40 

6 6.16 / 28.6 14.5 / 22.0 2.00 / 32.0 7.30 / 5.43 

7 6.46 / 28.7 12.5 / 21.1 3.40 / 32.5 5.90 / 4.70 

8 6.49 / 28.5 10.6 / 20.3 4.42 / 32.6 5.00 / 4.30 

TABLE I: Calculated values for the bulk modulus K and elas- 
tic constants of the BCC structure at several values of A = na, 
where a is the average interparticle spacing, for n — 0.23 (first 
number) and n = 0.50 (second number). All values are in 



units of N/r 



lattices. Their method relies on the fact that in inter- 
stitial regions, which occupy a large fraction of the sys- 
tem, the potential ^ is relatively uniform and changes 
abruptly in a boundary layer near the colloidal surfaces. 
As a result, the counterions are more tightly constrained 
to locations near the surfaces of these colloidal particles, 
and one can re-interpret these ions as part of an effec- 
tive sphere with a reduced charge. Again, since neither 
the effective charge nor the local ion concentration can 
be measured directly in experiments, the precise confir- 
mation of theory has proven difficult. Because our model 
allows the simultaneous prediction of the phase boundary 
and the shear moduli, we can avoid these uncertainties. 

There are, of course, additional interactions that we 
have neglected. The van der Waals attraction is much 
weaker than the screened Coulomb potential at the in- 
terparticle spacing of a sa 140 nm. We expect that cor- 
relation effects such as overcharging should have mini- 
mal consequence in a system of monovalent salt ions of 



NaCl . Dispersion forces are important at ion concen- 
trations greater than 100 /imol However, this effect 
occurs at distance scales on the order of 10 /im, and we 
do not expect this to change our model. 

Finally, we note that in the calculation of the bulk 
free energy, we employ free-volume theory even in the 
relatively low-particle-density regime - a regime in which 
pure hard-core interactions would predict a fluid phase. 
However, the screened Coulomb potential stabilizes the 
lattice structures and, in turn, lowers the melting density 
of the system. One could characterize the system with 
an effective density of higher magnitude, determined by 
the Barker-Henderson effective diameter of particles . 
A scheme like this would introduce another parameter 
into our model which would divide the screened-Coulomb 
potential into a "hard" part and a "soft" part. 

We have established a geometrical framework for un- 
derstanding the structural and mechanical properties of 
charged colloids. We have connected the phase diagram 
to the elastic moduli of the system with only one ad- 
justable parameter, the surface charge. The balance be- 
tween the drive of the system's entropy for maximal pack- 
ing fraction and the need of its screened repulsion for 
minimal interfacial area accounts for most of the prop- 
erties found in experiments. Its elucidation can serve as 
an intuitive guide to the engineering of all such charged 
systems. 
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